Magnetic neutron scattering from spherical nanoparticles with Néel surface anisotropy: analytical treatment

The magnetization profile and the ensuing magnetic neutron scattering signal from an inhomogeneously magnetized spherical nanoparticle with Néel surface anisotropy are derived analytically.


Introduction
Magnetic small-angle neutron scattering (SANS) is a powerful technique for investigating spin structures on the mesoscopic length scale ($1-100 nm) and inside the volume of magnetic materials Michels, 2021). Recent SANS studies of magnetic nanoparticles, in particular employing spin-polarized neutrons, demonstrate that their spin textures are highly complex and exhibit a variety of nonuniform, canted or core-shell-type configurations [see e.g.  Honecker et al. (2022) and references therein]. Surface anisotropy, vacancies or the presence of antiphase boundaries are generally considered to be at the origin of spin disorder in nanoparticles (Berger et al., 2008;Wetterskog et al., 2013;Nedelkoski et al., 2017;Kö hler et al., 2021a;Batlle et al., 2022). Magnetic SANS data analysis largely relies on structural formfactor models for the cross section, borrowed from nuclear SANS, which do not properly account for the existing spin inhomogeneity inside magnetic nanoparticles or nanomagnets (NMs).
Progress in magnetic SANS theory (Honecker & Michels, 2013;Michels et al., 2014;Mettus & Michels, 2015;Erokhin et al., 2015;Metlov & Michels, 2015Michels et al., 2016Michels et al., , 2019Mistonov et al., 2019;Zaporozhets et al., 2022) strongly suggests that, for the analysis of experimental magnetic SANS data, the spatial nanometre-scale variation of the orientation and magnitude of the magnetization vector field must be taken into account and macrospin-based models -assuming a uniform magnetization -are not adequate. The starting point for a proper analysis of the scattering problem is a micromagnetic continuum expression for the magnetic energy of the system. In the static case, this then leads to Brown's equations (Brown, 1963), a set of nonlinear partial differential equations for the magnetization along with complex boundary conditions on the surface of the magnet. From these equations the Fourier image and the magnetic SANS cross section may be obtained.
In this paper, we present an analytical treatment of the magnetic SANS cross section of a spherical NM with Né el surface anisotropy (Né el, 1954). This particular form of anisotropy arises because in an NM a significant fraction of atoms belong to the surface (with no neighbours on one side), and their magnetic properties such as exchange and anisotropy can be strongly modified relative to the bulk atoms.
The manuscript is organized as follows. In Section 2, we calculate the real-space spin structure of a spherical NM using classical micromagnetic theory within the second-order perturbation approach. In Section 3, we compute the threedimensional Fourier transform of the real-space spin structure, which directly yields the magnetic neutron scattering cross section and the pair-distance distribution function. The analytical results are benchmarked by comparing them with numerical finite difference simulations using the Landau-Lifshitz equation of motion. Finally, Section 4 summarizes the main findings of this study.
We also make reference to our accompanying numerical study (Adams et al., 2022) where, in contrast to the present analytical work, the full nonlinearity of the problem is considered.

Micromagnetic theory
In the static micromagnetic approach (Brown, 1963), the magnetic configuration of a system is described by the continuous magnetization vector field M(r), which has a constant magnitude kM(r)k = M 0 . The saturation magnetization M 0 is only a function of temperature. The normalized magnetization vector field is then defined as where r denotes the position vector. Our Hamiltonian for the NM includes the isotropic exchange interaction, the Zeeman energy, a uniaxial magnetic anisotropy for spins in the core and Né el surface anisotropy for those on the surface. In the continuum approach, it reads where A is the exchange stiffness constant, r is the del operator, Á is the Laplace operator, B 0 is a constant applied magnetic field, K c > 0 denotes the uniaxial core anisotropy constant, e A is a unit vector specifying the arbitrary core anisotropy axis and K s > 0 is the Né el surface anisotropy constant (Né el, 1954).
is the surface normal to the boundary of the NM, where and are the usual spherical angles (Garanin & Kachkachi, 2003;Kachkachi, 2007). In (2), the two surface integrals take into account the boundary conditions for the magnetization on the surface (@V) of the NM of volume V, which result from the exchange interaction and the Né el term. The magnetodipolar energy has been ignored in the calculations because of its mathematical complexity and since it is expected to be of minor relevance for smaller-sized NMs [see the recent atomistic simulations by Kö hler et al. (2021b)].
For small deviations from the homogeneous magnetization state, a perturbation approach is applicable. Let m 0 be the principal unit vector (average direction) associated with m(r) and let the vector function w(r) ? m 0 describe the spin misalignment. One can then write Assuming that x ; y ; z ( 1, the following second-order Maclaurin expansion in w is used to find an approximate closed-form solution for m(r): where m 0 is taken as a known constant vector in subsequent calculations. We choose the orthonormal vector base (Garanin & Kachkachi, 2009), and the parametrization and introduce the dimensionless coordinates n = r/R (with = knk = r/R), where r is the position vector, r ¼ ½r sin cos ; r sin sin ; r cos ; and R denotes the radius of the NM. The minimization of the Hamiltonian (2) then leads to the well known Helmholtz research papers equation with Neumann boundary conditions on the unit sphere (Kachkachi, 2007;Garanin & Kachkachi, 2003), where the constants are defined as with the dimensionless quantities The e (with = x, y, z) in (13) denote the unit vectors of the Cartesian laboratory coordinate frame (in which n and r are defined). We emphasize that there are only two independent differential equations for w, which is a consequence of the constraint km(r)k = 1. In our graphical representations, we will frequently use the following values: k c = 0.1 and k s = 3.0, which (using R = 5 nm and A = 10 À11 J m À1 ) correspond to K c = 80 kJ m À3 and K s = 12 mJ m À2 (Gradmann, 1986;O'Handley, 2000;Batlle et al., 2022). For M 0 = 1.7 Â 10 6 A m À1 , the relation between b 0 (dimensionless) and the external field is B 0 = (8/17)b 0 Â 1 T.
The fundamental solution of the homogeneous Helmholtz equation (9) is well known (Weber & Arfken, 2003;Riley et al., 2006). Its non-singular part can be expressed in spherical coordinates as an infinite series in terms of spherical harmonics Y 'm (, ) and spherical Bessel functions of the first kind j n (i ), The imaginary number 'i' in the argument of the spherical Bessel function is due to the negative sign in the Helmholtz equation (9). The expansion coefficients c 'm are obtained from the Neumann boundary condition (10) using the method of least squares (see Appendix A). From Appendix A it is seen that the zero-order term with ' = 0 vanishes. This physically makes sense, since the spin misalignment in our model is caused by the Né el surface anisotropy and thus, for symmetry reasons, there is no misalignment at the centre of the NM, i.e. ( = 0, , ) 0. By contrast, the largest spin misalignment is found at the boundary of the NM, i.e. = 1. Further, we find that the coefficients c 'm vanish in the case of odd ' and m, while they are real valued and even with respect to the index m, i.e. c 'm = c ';Àm . Taking these properties into account, one can conveniently express the solution in terms of the associated Legendre polynomials P m ' ðcos Þ with ' = 2 and m = 2 {note that we use the convention that Y 'm ð; with and In (18), , 0 is the Kronecker delta function, diag[ . . . ] denotes a 3 Â 3 diagonal matrix and Ç 0 ðÞ is the first-order derivative of (17) with respect to . For some small values of ' and m, the exact solutions of the integrals I 'm are listed in Table 1 in Appendix B.
From (18) it is seen that the functions depend linearly on k s , such that for k s = 0 the magnetization of the NM is homogeneous (as expected). Since we assume that x ; y ; z ( 1, it is clear that the validity of our solution is restricted to a finite range 0 k s k s, max . Taking only the terms with = 1 into account (corresponding to ' = 2), the remaining (second-order) expression reads ð; ; Þ ' À 15k s 32 where V V Vð; Þ ¼ diag cos 2 À 1=3 À sin 2 cosð2Þ cos 2 À 1=3 þ sin 2 cosð2Þ À2ðcos 2 À 1=3Þ 2 4 3 5 : A reasonable approximation for small in (21) is obtained by taking into account the first two terms in the infinite series (17) for Ç (). This results in the following expression [compare (21)]: In the limit !0 (small B 0 and small K c ), this expression reduces to a quadratic function in , The case of an infinite applied magnetic field B 0 , or of a strong uniaxial core anisotropy [compare (11) and (12) which recovers the expected result of zero spin misalignment. Note that the limit !1 is only obtained using all terms of the infinite series (17).
Of particular interest is the behaviour of as a function of the radius R of the NM. Inspecting the Hamiltonian (2), it becomes clear that the surface anisotropy energy scales as R 2 , while the uniaxial core anisotropy energy scales as R 3 . Since the core and surface anisotropies act in opposite ways (trying to make the spin structure more homogeneous and more inhomogeneous, respectively), we see that an increasing radius R corresponds to a decreasing . This behaviour reflects the NM's surface-area-to-volume ratio. With (21) it is not possible to make any prediction in this regard, because until this point we have not included the principal unit vector m 0 in the minimization of the Hamiltonian. Generally, m 0 is a function of k s , k c , b 0 and e A .

Figure 2
A comparison between the numerical solution using the Landau-Lifshitz equation (upper row) and the second-order analytical solution (21) for kwðnÞk ¼ ½ 2 1 ðnÞ þ 2 Fig. 1. This result was already predicted by Garanin & Kachkachi (2003). The solutions for 1, 2 (, , ) [using the particular m 0 (26)] then read 1 ' 15k s 32 2 ' 15k s 32 In Fig. 2, the analytical solution (21) (lower row) is compared with the numerical solution based on the Landau-Lifshitz equation (Bertotti, 1998), where is the gyromagnetic ratio, is the damping constant and the dot denotes the first-order time derivative [see our numerical study in the accompanying paper (Adams et al., 2022) for further details]. Shown is the vector norm of the w(n) function scaled to its maximum value. From Fig. 2 it is seen that our analytical approximation is in qualitative agreement with the results from the numerical simulation. The corresponding real-space spin structure m(n) is displayed in Fig. 3, where the surface spin disorder becomes clearly visible.
It is also instructive to compare our solution (21) with that obtained using the Green's function approach (Garanin & Kachkachi, 2003;Kachkachi, 2007). In particular, for located close to the surface, where the maximum spin misalignment with respect to m 0 occurs, the Green's function method yields the following approximate expression: This expression is also found when (21) is expanded in at the surface of the NM ( = 1). While the infinite series approach using spherical harmonics and spherical Bessel functions yields an exact solution of the Helmholtz equation, the Green's function approach provides an approximate explicit expression of in terms of the coefficients . Indeed, as was shown by Kachkachi (2007), in the presence of core anisotropy Green's function as the kernel of the Helmholtz equation is only obtained as a perturbative series in . As such, (29) is restricted to small values of , i.e. assuming that the core anisotropy and applied magnetic field are much smaller than the exchange coupling. This is manifest in (29) by the presence of the factor 1 À 2 =14 which implies that the contribution of spin misalignment may diverge when is too large (i.e. for a strong field and/or large core anisotropy).

Magnetic SANS cross section
The quantity of interest in experimental SANS studies is the elastic magnetic differential scattering cross section dAE M /d, which is usually recorded on a two-dimensional positionsensitive detector. For the most commonly used scattering geometry in magnetic SANS experiments, where the applied magnetic field B 0 k e z is perpendicular to the wavevector k 0 k e x of the incident neutrons (see Fig. 4 where V is the scattering volume and b H = 2.91 Â 10 8 A À1 m À1 is the magnetic scattering length in the small-angle regime (the atomic magnetic form factor is approximated by 1, The real-space spin structure in the x z plane computed using (4) and (21). Parameters are the same as in Fig. 2. The external field B 0 ' 266 mT is applied in the x z plane and inclined by an angle of = 40 relative to the z axis [compare with Garanin & Kachkachi (2003)].

Figure 4
A sketch of the perpendicular scattering geometry (B 0 ? k 0 ). The scattering vector q corresponds to the difference between the wavevectors of the incident (k 0 ) and scattered (k 1 ) neutrons. The angle q specifies the orientation of q on the detector. In the small-angle approximation, the component of q along k 0 is neglected. since we are dealing with forward scattering). e M MðqÞ = ½ e M M x ðqÞ; e M M y ðqÞ; e M M z ðqÞ represents the magnetization vector field M(r) in Fourier space, q denotes the angle between the scattering vector q and B 0 (not to be confused with the polar angle defined above), and the asterisk * stands for the complex conjugate. Note that in the perpendicular scattering geometry the Fourier components are evaluated in the plane q x = 0.
The Fourier transform of the three-dimensional magnetization vector field (with a tilde above the symbol) is defined as MðrÞ For subsequent calculations, we introduce the following dimensionless quantities: and we express the dimensionless scattering vector in spherical coordinates as Next, in (32) we use the following first-order approximation for the real-space magnetization vector m(n) [see (5) and (7)]: As shown in Appendix C, the final expression for the Fourier transform of the magnetization is then given by where and Ç ( ) is given by (17). The zero-order term / j 1 ()/ in (37) represents the form factor of a homogeneously magnetized sphere (Michels, 2021). In the limiting case of an infinite applied magnetic field, which is equivalent to the limit !1, the additional terms [second line in (37)] vanish [compare with (25)] and the spherical form factor remains. On the other hand, if k s = 0, the additional terms also vanish because, from the physical point of view, the Né el surface anisotropy cancels and from (18) we know that the coefficients a are linear in k s . Taking only the terms with = 1 into account and setting q = /2 ( x = 0), corresponding to the scattering geometry where the applied magnetic field B 0 k e z is perpendicular to the wavevector k 0 k e x of the incident neutrons (Fig. 4) where the radial function is R ðÞ can be approximated for small and, when only terms up to s = 1 in the infinite series (17) and (39) are kept, For small values, one finds the limit This can be seen by inspecting the definition of the Fourier transform in (32). Note that for q!0 the Fourier transform is proportional to the average of the magnetization vector field M and the maximum of this average is given by the homogeneous magnetization state. Using this result, the !0 limit for the first-order approximation in w of the Fourier transform of the magnetization yields Beyond the linear approximation in w, a non-vanishing term appears in f M M M M M M in the limit !0, which reduces the Fourier components relative to the homogeneous magnetization state. In the second order in w, the result is [compare (5)] Using (34) and the dimensionless two-dimensional magnetic SANS cross section S M ð; q Þ can be straightforwardly obtained as [compare (31)] In the limit k s !0, the resulting cross section from (37)  m 2 0;x þ m 2 0;y cos 2 q þ m 2 0;z sin 2 q À À 2m 0;y m 0;z sin q cos q Á : Relation (49) nicely demonstrates that, depending on the orientation of the uniformly magnetized particle, different angular anisotropies become visible on the detector. For m 0 k e x (i.e. m 0y = m 0z = 0) the scattering pattern is isotropic, while it exhibits a cos 2 q (sin 2 q ) type shape when m 0 k e y (m 0 k e z ). Fig. 5 shows S M ð; q Þ along with the contribution of the individual Fourier components to (48). The upper row in Fig. 5 presents the results taking into account only the zeroorder term [j 1 ()/]m 0 from (40), while in the lower row the second-order term ( = 1) is additionally included. Since the zero-order term represents the case of a homogeneously magnetized NM, this comparison provides useful insights about the impact of the Né el surface anisotropy on the magnetic SANS cross section. In the case of a uniformly magnetized NM (upper row) the Fourier components j f M M x j 2 , j f M M y j 2 and j f M M z j 2 are isotropic (rotational symmetry), while including the second-order terms (lower row) leads to anisotropic behaviour of the transverse components j f M M x j 2 and j f M M y j 2 . The cross term (CT) averages to zero for both situations and the dominant contribution to the magnetic SANS cross section (for the parameters chosen in Fig. 5) is given by the j f M M z j 2 component. Therefore, it may be concluded that the impact of the Né el surface anisotropy on S M ð; q Þ is relatively small. By comparing the S M ð; q Þ from the upper and lower rows, it is seen that by including the Né el surface anisotropy the circular symmetry of the zeros of S M (deep-blue colours) is broken. This feature becomes more clearly visible by analyzing the azimuthal average of S M ð; q Þ, which is readily computed as In the limit k s !0, the azimuthal average corresponding to (49) is lim k s !0 We have also calculated the pair-distance distribution function and the correlation function In the limit k s !0, the pair-distance distribution and the correlation function corresponding to (51) are lim k s !0 These functions are displayed graphically in Fig. 6. Due to the surface-anisotropy-induced spin disorder, the form-factor extrema of IðÞ [ Fig. 6(a)] are damped and shifted slightly to larger q values [i.e. smaller structures; compare the first minimum of IðÞ for k s = 3]. Moreover, as already observed in numerical micromagnetic continuum simulations (Vivas et al., 2017(Vivas et al., , 2020, the oscillations are damped for the case of surface spin disorder, which mimics the effect of a particle-size distribution or of instrumental resolution. In agreement with  (40), which corresponds to the case of a homogeneously magnetized particle. The lower row displays the results when the second-order term ( = 1) in (40) is taken into account. The parameters are e A = e z , b 0 = 0.1e z (B 0 ' 48 mT), k c = 0.1, k s = 3 and m 0 = [sin cos , sin sin , cos ]. Note that y and z denote the dimensionless components of the scattering vector [compare equation (34)]. Since the Né el surface anisotropy effectively has a cubic symmetry (see Fig. 1), we average S M over the angles = (45 , 135 , 225 , 315 ) and = 20 . A logarithmic colour scale is used. this observation is the finding that the maximum of the PðÞ function [ Fig. 6(b)] appears at smaller distances than in the homogeneous case. Likewise, due to spin disorder, the CðÞ function [ Fig. 6(c)] exhibits a larger amplitude (Mettus & Michels, 2015).
To analyze the role of the surface anisotropy more quantitatively, we have computed the following quantities, which describe the deviation of the one-dimensional SANS cross section and of the pair-distance distribution function from the homogeneous particle case: Pðk s Þ ¼ Fig. 6(d) depicts both Iðk s Þ and Pðk s Þ as a function of k s . The difference is only of the order of a few percent, which suggests that the effect of surface anisotropy on the SANS observables is relatively weak within the present analytical approximation: see our accompanying numerical work (Adams et al., 2022), which takes into account the full nonlinearity of the micromagnetic equations. However, this is only true for the magnetic interactions considered here. Taking into account the anisotropic and long-range dipole-dipole interaction and the asymmetric Dzyaloshinskii-Moriya interaction will very likely result in more inhomogeneous spin structures and in larger deviations from the macrospin model (Vivas et al., 2017(Vivas et al., , 2020Pathak & Hertel, 2021). Likewise, for NMs of elongated shapes, the surface anisotropy renders an additional first-order contribution to the effective energy (Garanin & Kachkachi, 2003), in addition to the second-order cubic contribution discussed above. This new shape-induced contribution could also lead to an enhancement of the spin misalignment. The analytical calculations presented here provide a general framework for future studies of more complicated (anisotropic) magnetic interactions. The approach can be straightforwardly adapted to other particle shapes such as a circular planar disc.

Conclusions
We have analytically computed the magnetization distribution and the ensuing magnetic SANS cross section of a spherical nanoparticle. Our micromagnetic Hamiltonian takes into account the isotropic exchange interaction, an external magnetic field, a uniaxial anisotropy for the particle's core and Né el anisotropy on its boundary. The resulting Helmholtz equation has been solved by expanding the real-space magnetization in terms of spherical Bessel functions and spherical harmonics. The central results are the infinite series (16) and its second-order expansion (21) for the real-space magnetization, and the corresponding Fourier transforms (37) and (40). Using these expressions, the two-dimensional magnetic SANS cross section S M ð; q Þ, the azimuthally averaged SANS signal IðÞ, and the correlation functions PðÞ and CðÞ have been obtained and compared with the case of a homogeneous spin configuration (uniform magnetization vector field). The signature of Né el surface anisotropy (of constant k s ) has been identified in all of these functions. However, its effect is relatively small, even for large values of k s . Taking into account the magnetodipolar and/or the Dzyaloshinskii-Moriya interaction, or shape asymmetry, will probably result in configurations with stronger spin misalignment (e.g. in vortex-type textures or skyrmions) and thereby in more prominent signatures in the SANS cross section and correlation function. These interactions are beyond the scope of the current analytical approach and will be considered in our future (numerical) work.

APPENDIX A Solution of the boundary value problem of the Helmholtz equation
The coefficients c 'm in the fundamental solution (15) of the Helmholtz equation (9) must be determined such that the Neumann boundary condition (10) is satisfied. For this purpose, we use the method of least squares, where we make use of the orthogonality properties of the spherical harmonics Y 'm (, ). The normal derivative of (15) where Our goal is now to minimize the following error functional with respect to the coefficients c 'm : The minimum of this error functional is found from the condition that the partial derivatives of "½c 'm with respect to ðc ij Þ Ã vanish, @"½c 'm @ðc ij Þ Ã ¼ 0; where * denotes the complex conjugate. Using the orthogonality relation (62) (63), the solution of (61) can be written as Alternatively, again using the indices ' and m, expressing the coefficient as in (13) and using the matrix vector product, the coefficients c 'm can be more explicitly written as For some low orders of ' and m, the exact solutions of the integrals I 'm are presented in Appendix B. From (65) several conclusions can be drawn. First, the zeroorder term (' = m = 0) in (15) vanishes, which can easily be shown by rewriting the coefficient c 00 as and this is due to the orthogonality of m 0 and g with 2 {1, 2} [see (6)]. Moreover, we see that which is a consequence of the behaviour of the spherical Bessel functions of the first kind at the origin. Since the j ' with ' ! 1 all vanish at the origin = 0, this implies that the total also vanishes at the origin. Note that j 0 (0) = 1 but does not contribute to since c 00 = 0. From the physical point of view this makes sense, since the spin misalignment is caused by the Né el surface anisotropy and thus, for symmetry reasons, there is no spin disorder at the centre of the spherical NM; the highest misalignment is found at its surface. Secondly, in Table 1 in Appendix B it is seen that the coefficients I ';m vanish for odd ' or m and that I ';m = I ';Àm , so that the expansion coefficients also exhibit this symmetry, Taking these properties into account, one can express the solution (15) more conveniently in terms of the associated Legendre polynomials P m ' ðcos Þ with ' = 2 and m = 2 {note that we use the convention that Y 'm ð; [p. 378 (14.30.1) of Olver et al. (2010)]}, where Ç () is defined in (17) and a in (18). The infinite series in (18) The integrals U ð1Þ m , U ð2Þ m and U ð3Þ m are solvable straightforwardly using Euler's formula for the complex exponential and by splitting the region of integration according to the absolute values of the trigonometric functions. In the denominator of U ð1Þ m and U ð2Þ m we have included the Kronecker delta |m|, 1 to take account of the cases m = AE1. It is common to express the integrals K ð1Þ 'm and K ð2Þ 'm by the substitution x = cos (dx = Àsin d), Using U ð3Þ m ¼ 2 0;m , we need only compute K ð2Þ ';0 such that the associated Legendre polynomials in K ð2Þ ';m are reduced to the Legendre polynomials with one index only [p. 352 of Olver et al. (2010)]. By considering the symmetry properties of the Legendre polynomials it becomes clear that the integral must vanish for odd ' and can be simplified for even ' in the following way: The closed form of K ð2Þ ';0 is then found in terms of the Gamma function [p. 771 (7.126.1) of Gradshteyn & Ryzhik (2007) The overall solution for I z 'm is then written as  Table 1 Values of the integrals (19) for some small values of ' and m.
From these results it is seen that the integrals I 'm with 2 {x, y, z} vanish for odd ' and m [note that in (86) this is only the case for m = 2]. For the remaining integrals K ð1Þ 2;2 in (86), where ' = 2, we do not give an expression in closed form. However, one must exist in terms of the Gamma function or the Beta function, since [p. 324 (3.251.2) of Gradshteyn & Ryzhik (2007)] where B(Á, Á) is the Beta function (Euler integral), and the associated Legendre functions P 2 2 of even order and degree are true polynomials, as seen for example from the related Rodrigues formula [p. 360 (14.7.14) of Olver et al. (2010)]. We used Mathematica (Wolfram Inc.) to determine the integrals up to the sixth order in ' (see Table 1).

APPENDIX C Derivation of the Fourier transform of the magnetization
The Fourier transform of the magnetization vector field M(r) is written as In the following, we will use dimensionless quantities. For this purpose, we define the dimensionless scattering vector t = qR, where R is the radius of the nanomagnet, the dimensionless position vector n = r/R and the dimensionless magnetization vector m = M/M 0 , where M 0 is the saturation magnetization. Substituting in (88) The next step consists of calculating the Fourier integral of the first-order approximation (36) of the magnetization vector m.
Since (36) is formulated in dimensionless spherical coordinates , , , it is convenient to express the scattering vector t in spherical coordinates as well, n ¼ sin cos ; sin sin ; cos ½ ; ð92Þ t ¼ sin q cos q ; sin q sin q ; cos q Â Ã ; so that the plane-wave expansion of the complex exponential can be used (Jackson, 1999), exp ðÀit Á nÞ ¼ 4 P 1 k¼0 P k n¼Àk ðÀiÞ k j k ðÞ Y Ã kn ð; Þ Y kn ð q ; q Þ: The Fourier integral (90) We now use the infinite series (15) for to express the firstorder approximation of the magnetization, which leads to Since the integral transform (96) is linear, each term in (97) can be separately transformed. For the zero-order term, we obtain This result is well known to the neutron-scattering community as the spherical form factor, corresponding to a uniformly magnetized spherical particle (Michels, 2021). In the second step, we carry out the integration of the higher-order terms from (97). The radial and angular parts in the higher-order terms of (97) are multiplicative, such that the volume integral (96) is separable into As The integral (100) directly corresponds to the orthogonality relation of the spherical harmonics and thereby we have B kn 'm = 'k mn [p. 378 (14.30.8) of Olver et al. (2010)]. Due to the term 'k in B kn 'm and since B kn 'm and A 'k are multiplicative in (101), the following spherical Hankel transform results: